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ABSTRACT 

We analyse a 10,000 star spectroscopic survey, focused on Galactic thick disk 
stars typically 2-5 kpc from the Sun, carried out using the AAOmega Spectrograph 
on the AAT. We develop methods for completeness-correction of the survey based on 
SDSS photometry, and we derive star distances using an improved isochrone-fitting 
method with accuracies better than 10%. We determine the large-scale kinematic 
(V0=172 km/s and a{4>,r,z) = (49,51,40) km/s.), and abundance properties of the 
thick disk, showing these representative values are a fair description within about 
3 kpc of the sun, and in the range 1-3 kpc from the Galactic Plane. We identify a 
substantial overdensity in lines of sight towards the inner Galaxy, with metallicity 
[Fe/H]'^ —1 dex, and higher line of sight velocities than the thick disk, localised 
along the direction {£,b) = (48°, —26°). This overdensity appears to be towards, but 
closer than, the known Hercules-Aquila halo overdensity and may be related to the 
Humphreys-Larscn inner galaxy thick disk asymmetry. 
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1 INTRODUCTION 

The advent of large-area photometric and spectroscopic sur- 
veys in the recent decade has made observational studies of 
the Galactic stellar populations extremely vibrant. Rapid 
progress is being made in obtaining independent descrip- 
tions of the chemical, kinematic and spatial distributions, 
especially of the thick disk and halo, although systematic 
challenges remain associated particularly with the defini- 
tions of stellar populations. Photometric definitions, which 
led to discovery of th e thick disk, are reviewed by, for ex- 
ample, ICilmore et al.l ([l98S ) and the history of suc h stud- 
ies up to the present is described by lYoshiil (|2013l ) . Much 
recent focus has been on chemical abundance definitions. 
Those studies based on either small numbers of high quality 
spectra or large numbers of lower resolution spectra, have 
become a major industry. They have established that the 
thick disk is a discrete stell ar populatio n , at l east in the 
solar neig hbourhood (but cf. IBovv et al.l (|2012h for an or- 
thogonal view), but raised many open questions and incon- 
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sisten t resul t s on gr adients, origins, a nd so on JPuhrmannl 
(|2004 l2008l. l201ll). iFuhrmann et al.) (l2012t). iReddvl (1201 
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2007 T.lReddv fc LamberTi 200^ . iReddv et al I (1200a 
Ben sbv fc Feltzind ('2012a tf ). iBensbv et al.l (l201ll. 
Lee et all (120111). [Schlesingeret all (1201 ll). ICheng et al 
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As just two mor e of these very many recent examples_ 
lJuric et al.l (|2008l ) determined the exponential scale height 
and scale length of the thick disk from Sloan Digital Sky 
Survey (SDSS) photometric d a ta to be 0.9 kpc and 3.6 
kpc respectively. ICarollo et al.l (|201(]| 'l re-visited the kine- 
matics of the thick disk and halo using SDSS Data Release 
7 (DR7) and particularly investigated the overlaps between 
known stellar populations to conclude that the metal weak 
thick disk is a kinematically and chemically distinct compo- 
nent between the halo and thick disk. The main challenge 
for these kinds of studies has been the intrinsic overlap of 
population distribution functions, and the associated com- 
plexities in how to deconvolve them correctly. Improving 
constraints then requires a moderately high surface density 
spectroscopic survey, with accurate photometric and proper 
motion data, sampling a sufficient range of Galactic coor- 
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dinates to support statistical separation of the populations. 
Much is being learnt, but much remains to be earned, es- 
pecially from spectrosc opy outside the Solar neighbourhood 
(for example, Kordop atis et. al.l l|201ll )') 

We present here analysis of a 10,000-star spectroscopic 
survey, with Anglo- Australian Observatory AAOmega spec- 
troscopy complementing SDSS photometry, that has the 
potential for single-population and correlated population- 
overlap studies of the kinematics, chemistry and spatial dis- 
tribution of the thick disk. This survey uses F/G stars to 
probe significant distances from the Sun; the chemical abun- 
dances of these stars provide information on the enrichment 
history during early stages of galaxy formation. The radial 
velocities, along with chemical abundance and spatial distri- 
bution, allow us to differentiate between populations. The 
main scientific objective of this survey was to determine the 
unbiased distribution functions over kinematics, metallicity 
and spatial distributions of the thick disk and halo, with 
emphasis on the region of overlap between these two popu- 
lations. 

In section 2 of this paper we present the data, dis- 
cuss the selection function and compare the photometry 
and spectroscopy to SDSS DR7. In section 3 we present 
a weighting method for completeness correction, and there- 
after explain an isochrone fitting technique used to derive 
distances to the stars. In this same section, we also pro- 
vide a brief discussion of the proper motions (taken from 
SDSS DR7) and show that the resulting tangential veloci- 
ties have inevitable large uncertainties. We hence show that, 
for this study, it is better to exclude proper motions from 
the analysis and model the data with parameters derivable 
from the spectroscopic measurements i.e. line of sight ve- 
locities, metallicities and distances. In section 4 we outline 
our method to implement such modelling by marginalising 
over unknown parameters. We then derive and present our 
results for the best-fitting kinematic and photometric model 
of the thick disk, and analyse and characterise the residuals 
of this model when applied to the data, searching for - and 
identifying - significant phase-space substructure. 
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Figure 1. The WGN sel ection function to isolate stars from the 
thick disk-halo interface l|Wvse et al.ll2013h [ al colour-magnitude 
diagram of data from SDSS DR7 covering 32 fields spanning 
-1° < Dec < 1° and 133° < RA < 339°, and a broad selection 
function ; 14 ^ go < 22, -1 (3 - r)o ^ 1.7, -1 s; (r - i)o sS 1, 
(u — g)o ^ 0. The main Galactic stellar populations seen in this 
CMD are - the thin disk near {g — r)o ~ 1.3, the thick disk 
(brighter) and halo (fainter) near {g — r)o ~ 0.4; [b] CMD of (a) 
with the following selection : —0.5 < (r — i)o < 0.6, {u — g)o ^ 0.6. 
This significantly reduces thin disk contamination. The area en- 
closed within the dot-dashed box denotes the region that is se- 
lected from g vs {g — r)o space to isolate stars from the thick 
disk-halo interface; [c,d] Colour-colour diagram for the same data 
as in (a). The {u — g)o cut mainly eliminates quasars while the 
(9 ~ '')o and (r — i)o cuts eliminate thin disk and binary stars. 
The dot-dashed box again indicates the selected region of thick 
disk-halo interface stars isolated for study here. 



2 THE DATA 
2.1 The Survey 

A wide-fiel d spec tr oscop ic survey was carried out by 
IWvse et al.l (l2006l . 12011 (hereafter WGN)) using the 
AAOmega Spectrograph on the Anglo Australian Telescope, 
fed by the 400-fibre 2dF facility, targeting lines of sight 
in the thick disk-halo interface, each covering a two de- 
gree field of view. Spectroscopy was obtained for ~ 10, 000 
stars, drawn from the Sloan Digital Sky Survey (SDSS) pho- 
tometric catalogue, along 32 lines of sight in the region 
133° < RA < 339°; -1° < Dec < 1° i.e. along the celes- 
tial equator. These cover Galactic longitudes between 228° 
and 67° at intermediate latitudes (26° < |6j < 63°). Table 
[1] provides the coordinates of all the observed lines of sight 
along with the number of stars analysed in each field. 

Iron abundances were derived from medium resolu- 
tion spectra (R=6000) covering 3700A to 4700A us- 
ing a meth o dolog y based originally on that developed by 
iBeers et aU (| 19991 ) and extended by WGN. The technique. 



which uses the CaK-line to determine iron abundance, is 
calibrated using halo and thick disk field stars. We mea- 
sured abundances only for spectra having more than 100 
net counts per 0.34A pixel at 4150A. The typical median 
counts for the sample is 250 per 0.34A pixel at 4150A. 
As noted below, this leads to an abundance uncertainty of 
A[Fe/H] ~ 0.2 dex. Radial velocities were determined by 
cross correlation against sky and standard stars, providing 
velocities with random plus systematic accuracy better than 
10 km/s. 

The selection function employed in obtaining the data 
was as follows: 

17.0 <g < 18.5 
-1.0 < {g-r)o < 0.5 
-0.5 < (r - i)o < 0.6 
(u - g)o > 0.6 

This selection function preferentially selects thick disk and 
halo stars, mainly F/G dwarfs, and minimises contamination 
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Table 1. Fields Covered by the WGN Data 



RA,Dec 


l,b 


Mean Reddening 


Number of Stars 


Number of Stars 


% Completeness of 


deg 


de 


g 


E(B-V) 


analysed from WGN 


SDSS DR7 


WGN spectra wrt SDSS DR7 


133 


228, 


27 


0.033 


288 


2907 


10 


135 


229, 


28 


0.035 


236 


2523 


9 


139 


231, 


32 


0.030 


254 


1855 


13 


141 


233, 


33 


0.032 


26 


1713 


\ 


145 


235, 


37 


0.064 


268 


1524 


17 


153 


242, 


43 


0.036 


255 


1280 


20 


159 


247, 


48 


0.065 


304 


1210 


25 


161 


249, 


49 


0.047 


299 


1238 


24 


163 


251, 


50 


0.046 


130 


1211 


10 


175 


268, 


58 


0.021 


134 


1187 


11 


179 


275, 


60 


0.022 


118 


1255 


9 


185 


286, 


62 


0.024 


262 


1242 


21 


189 


295, 


63 


0.021 


123 


1376 


9 


197 


312, 


63 


0.023 


294 


1559 


19 


199 


316, 


62 


0.028 


253 


1546 


16 


203 


324, 


61 


0.027 


274 


1717 


ifi 
±\j 


205 


328, 


60 


0.028 


257 


1790 


14 


213 


341, 


57 


0.041 


199 


2190 


9 


217 


347, 


54 


0.039 


270 


2428 




223 ,0 


354, 


50 


0.047 


288 


2836 


10 


229 ,0 


1, 45 


0.057 


574 


2669 


15 


235 ,0 


6, 41 


0.100 


209 


4668 


5 


237 ,0 


7, 39 


0.093 


275 


4859 


5 


239 ,0 


9, 38 


0.113 


447 


4401 


10 


249 ,0 


16, 29 


0.107 


600 


6150 


10 


313 ,0 


48, 


-26 


0.096 


756 


10754 


7 


315 ,0 


49, 


-28 


0.081 


498 


9101 


5 


317 ,0 


50, 


-30 


0.097 


235 


7614 


3 


319 ,0 


51, 


-31 


0.081 


527 


6410 


8 


323 ,0 


54, 


-35 


0.048 


194 


4732 


4 


333 ,0 


62, 


-43 


0.069 


264 


2477 


10 


339 ,0 


67, 


-47 


0.063 


273 


1680 


16 



from the thin disk population. In addition, in several fields, 
fainter halo turnoff candidate stars were targeted to exploit 
periods of very good seeing, together with Blue Horizontal 
Branch stars. This extra selection was 

18.5 < g < 19.5 
{g - r)o < 0.25 

although in practise few stars in this faint extension 
delivered spectra of sufficient signal-to-noise ratio to deliver 
abundances, and essentially none with g > 19.0. 

Figure [TJa) shows the colour-magnitude diagram 
(CMD) for all stars in the region 133° < RA < 339° 
and -1° < Dec < 1° from the SDSS Data Release 7 
l|Abazaiian et al]|2009l ). The group of stars centred around 
{g—r)o = 1.3 and faint apparent magnitudes are intrinsically 
faint, nearby, low mass stars. These, from their kinematics 
and metallicity estimates, are identified to be thin disk stars. 
The stars centred around {g — r)o = 0.4 are older popula- 
tions of stars with a bluer main sequence turn-off. There 
is a shift at about go = 18.5 in the average colour of the 
stars, and this is caused due to the different metallicities of 
the two groups of stars. The stars with brighter apparent 
magnitudes have a relatively redder turn off, due to rela- 
tively higher metallicity - this is the thick disk population. 
The stars with bluer turn-off have comparably old ages, but 



lower metallicity compared to the thick disk, and are halo 
stars. 

Figure [TJb) illustrates how the CMD changes when the 
(r — i)o and (u — g)o colour cuts are imposed on[TJa). The 
number of stars is reduced to 75% of the original sample. The 
(r — i)o colour-cut preferentially picks out the thick disk and 
halo stars, over the thin disk ones, while the {u — g)o colour 
cut eliminates most quasars. This fact is well illustrated by 
the colour-colour diagrams in Figure[TJc) and (d). The black 
squares in[TJb) mark out the g and {g — r)o colour cuts im- 
posed in this WGN AAOmega Survey. This clearly shows 
that the selection function achieved the goal of preferentially 
identifying thick disk and halo stars, while minimising con- 
tamination from the thin disk and quasars. The F/G dwarfs 
thus selected are at distances of a few kpc. 

Data Release 4 (DR4) of SDSS served as the original in- 
put catalog to the AAOmega Survey. Thus the astrometric 
and photometric information of the WGN selection func- 
tion is from SDSS DR4. For this work, the photometry was 
updated to the more recent Data Release 7 (SDSS DR7). 

The stars were cross-matched between DR4 and DR7 
based on their unique SDSS identifier. Astrometry and pho- 
tometry compare well between SDSS DR7 and DR4. The 
positions of stars from the two datasets are consistent, as 
expected. The extinction corrected apparent magnitudes 
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Figure 2. Comparison of derived quantities from Spectroscopy 
from WGN (AAOmcga Survey) to the cross-matched spec- 
troscopy from SDSS DR7. 84 matches are found based on the 
position of stars. The top 2 panels show scatter plots of [Fe/H] 
and radial velocity from WGN against the values from SDSS 
DR7, where the dashed line indicates the one-to-one correlation 
for WGN data. The bottom panels show the distribution of the 
difference between the two measurements. The black histogram 
denotes the difference between WGN and SDSS DR7, the blue 
smooth line shows a gaussian fit to this histogram. The dashed 
vertical line in each case marks the mean offset of SDSS DR7 
spectroscopy from WGN, while the dotted vertical lines marks 
the 1(7 limits of the scatter. The value of the mean offset and 
1(T limits for both [Fe/tl] and radial velocity is indicated in the 
respective plots. 



from DR4 are displaced from the corresponding values from 
SDSS DR7 by 0.002 magnitude with a la spread of 0.01 
magnitude, which is negligible. Similarly so for the g — r 
colour. This suggests that the recalibration of photome- 
try/extinction effects in DR7 as compared to DR4 is mini- 
mal. 

Figure [2] shows the comparison between the spectro- 
scopic results of WGN (obtained from measurements in the 
AAOmega survey) for the 84 stars in common with the spec- 
troscopic observations released as part of SDSS DR7. Metal- 
licity values between the two surveys are in very good agree- 
ment: they are offset in the mean by less than about 0.1 dex 
and show a scatter of 0.2 dex at the la level. Radial veloci- 
ties from the two surveys are also in good agreement, with 
a mean offset less than 2 km/s (la = 3 km/s). (The few 
large outliers are likely to be due to variable stars). Both 
the metallicity and velocity offsets indicate uncertainties we 
will adopt for the AAOmega survey. 

We thus have available for analysis 32 fields spanning 
the equatorial stripe, each some 2 degree squared area, cen- 
tred on the listed line of sight coordinates, with [Fe/H] and 
radial velocities from WGN and ugriz photometry and po- 
sitions from SDSS DR7. 

Figure |3] presents an overview of the spectroscopic mea- 



surements from the WGN survey. It shows a plot of the 
measured heliocentric velocity and metallicity vs RA. This 
broadly indicates the presence of both the non-rotating, 
metal-poor halo component and the rotating and relatively 
metal-rich thick disk component. An azimuthal velocity of 
180 km/s at a distance of 2 kpc projected onto the line of 
sight is over-plotted for reference. 



3 PROPERTIES AND CALIBRATION OF THE 
WGN SAMPLE 

3.1 Completeness Correction 

Completeness is always a significant issue with any spec- 
troscopic dataset as there are many more stars which are 
present in a given direction of the sky than the number for 
which we have spectra. With the introduction of multi-fibre 
spectrographs in the last decade, spectra of multiple objects 
are now obtained simultaneously but this is still restricted 
to a few hundred objects per field of view. 

There is a further bias towards brighter stars which 
arises from discarding those spectra from the sample which 
have low signal-to-noise ratio. The sample for which spec- 
trometric parameters are derivable has an additional bias 
against metal-poor, hot, blue stars, which generates biased 
completeness in the purely photometric input selection func- 
tion, described in Section 2.1. 



3.1.1 Completeness of SDSS DR 7 

We first check the completeness of SDSS by comparing to a 
deeper dataset, namely COMBO-17. COMBO-17 is a spec- 
trophotometric survey which imaged 5 fields at high Galactic 
latitudes covering 1 square degree of the sky in 17 optical 
filters. It obtained deep and sharp r-band data and the sele c- 
tion of stars is complete down to r ~ 23 (|Wolf et al.ll2004l ) 
The SDSS DR7 photometry is sufficiently close to 100% 
complete in the magnitude range 16.0 < r < 21.0 (which 
encompasses the range of the WGN dataset 16.5 < r < 
19.0), in comparison to overlapping fields in the COMBO- 
17 survey. Q. 

It is thus reasonable to adopt SDSS DR7 as complete 
in the magnitude range of relevance to WGN. We therefore 
go on to use SDSS DR7 to evaluate the completeness of the 
WGN sample. 



3.1.2 Completeness of WGN 

We find that the WGN dataset has completeness levels that 
are as high as 25% of all potential candidates having spec- 
tra of sufficient quality to allow derivation of both [Fe/H] 
and radial velocity. The last 2 columns of Table [T] show the 
number of stars in the given field in the SDSS master pho- 
tometric catalog (employing the same selection functions as 
in the data) and the relative completeness of the WGN data 
field- by- field. 
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the comparison of complete- 

between SDSS and COMBO-17 see 
http: / / www.sdss3.org/dr8 / imaging/ other _info.php 
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Figure 3. Measured heliocentric velocity and metallicity as a function of RA from the WGN dataset analysed here. The corresponding 
Galactic longitude / and latitude b of the fields is shown in parallel axes. RA, I and b are quoted in degrees. The colours indicate the 
number of stars within each pixel (colour-coded in the colour bar) going from red for small to blue for large star count. The black dashed 
line shows an azimuthal velocity of 180 km/s at a distance of 2 kpc projected onto the line of sight and is over-plotted to guide the eye 
along the thick disk population. 



3.1.3 Weighting Method 

We describe here a weighting method in colour-magnitude 
space to correct the WGN spectroscopic sample for com- 
pleteness inside the WGN selection function. This method 
is widely applied in sampling, in situations as diverse as 
opinion polling and ornithology. We note its applications in 
studies of the thick disk in astronomy below. 

Stars in a colour-magnitude bin of the spectroscopic 
data are weighted up to the actual number of stars available 
for observation in that bin, by a weighting factor W which is 
obtained by taking the ratio of the density function of stars 
in colour magnitude space D[go, (g — r)o] of the SDSS DR7 
catalogue to the data themselves: 



W = 



Ds 



7 [go, (g - l-)o] 



DwGN[go,(g-r)o] 



(1) 



Three different methods of computing the density func- 
tion for both the Sloan photometric superset and the WGN 
dataset were tried : (a) probability density functions were 
calculated by convolving each data point with a Gaussian 
kernel function; (b) one dimensional number density distri- 
butions (as a function of apparent magnitude) for individ- 
ual colour bins were determined by simple polynomial fits to 
ID apparent magnitude histograms; (c) a two-dimensional 
histogram (i.e. number of stars as a function of colour- 
magnitude) was used to represent the distribution of stars 
in the colour-magnitude plane. 

The three methods provide progressively decreasing de- 



grees of smoothing ((a)-most smoothing, (c)-least smooth- 
ing) across the colour-magnitude plane. 

The pixel or bin choices in each of these methods was 
made so that we sufficiently resolve the main features of 
interest on the CMD (i.e. the different Galactic stellar pop- 
ulations), while ensuring that the pixels contain sufficient 
number of stars relative to the noise level. This choice may 
also be dependent on line-of-sight. 

We found that the distribution functions produced by 
these three methods did not differ significantly, indicating 
that the particular method chosen to compute the density 
distributions in colour magnitude space is not crucial. We 
chose method (c) as the most robust of the three methods 
for the rest of the analysis. Along with best performance 
with low number statistics and ease of computation, this 
method has particular benefits over the other two. Method 
(a) involves the risk of the smearing of tiny features that may 
be of importance through over-smoothing, while edge effects 
introduce problems in the polynomial fitting technique of 
method (b). 

This weighting function was computed field by field, 
taking into account that extinction varies across the sky 
and that the form of the CMD also varies as an effect of 
the different density profiles of the halo and thick disk. The 
complete photometric sample from SDSS DR7 was selected 
for the same colour magnitude cuts, but for a larger spa- 
tial region (6x6 degrees) centred on each field, to minimise 
noise. Colour-magnitude bins of size Ag, A{g — r) = 0.7, 0.13 
mag were found to be an appropriate choice for all fields. 
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Figure 4. Completeness correction is implemented for the WGN data through a weighting method. This method is illustrated here for the 
(not untypical) field with coordinates RA,Dec = 203°, 0°. Left: the distribution in colour magnitude space of the stars observed by WGN in 
the 2x2 degrees field. Middle: SDSS DR7 data for a 6x6 degrees region centred around the same coordinates, for the same colour-magnitude 
cuts. Right: the corresponding weights, which are calculated as the ratio of the middle to the left-most panel, multiplied by the area factor, 
which is given by the area covered on the sky by the left-most panel divided by the area covered by the middle panel 



The weighting method is illustrated in Figure 3] for a sin- 
gle field, centred at RA,Dec= 203°, 0°. The first panel shows 
the WGN data for a single 2dF field (i.e. a circular field of 
radius 1 degree centred on the mentioned coordinates) . The 
middle panel shows the SDSS DR7 data for a larger region 
centred around the same coordinates, for the same colour- 
magnitude cuts. The last panel shows the weights, which are 
calculated as the ratio of the counts in the middle to those 
in the left-most panel, multiplied by an area factor (given by 
the area covered on the sky by the left-most panel divided 
by the area covered by the middle panel). In order to avoid 
large errors from low signal-to-noise pixels, the distribution 
of weights in every field is cut off at the 90th percentile. Ob- 
jects having weights greater than the 90th percentile value 
are reassigned to this new maximum value. 

In Figure[S]we illustrate how the metallicity and velocity 
probability distribution functions (PDF) change when the 
completeness correction is applied. The histograms shown 
are for the high-latitude fields, with RA from 190° to 220°. 
This covers 6 fields at high latitudes (54° < b < 63°) in the 
inner galaxy (312° < I < 347°). The data have been split 
into faint {go > 18, top panels) and bright {go < 18, bottom 
panels) stars. In both cases, the unweighted metallicity dis- 
tributions show a high number of metal-poor stars, particu- 
larly in the faint bins, which are likely to be predominantly 
halo stars given these are high latitude fields. Since we are 
looking towards the inner galaxy, metal-rich thick disk stars 
must also be numerous. Applying the weights corrects the 
distribution for the fact that the data do not show high 
enough normalisation at the metal-rich end, the effect being 
more evident in the faint bin. We note that the heliocentric 
velocity distribution is centred on zero as one would expect 
for the chosen line-of-sight towards the inner Galaxy. We can 
see that the velocity distributions are in general smoothed 
by the completeness correction. 

3.1.4 Errors on the Weights 
The error on the weight W was calculated as 

W 

a» [go, (g - r)o] = , ^ ^ (2) 

VDwGNlgo, (g - r)oJ 



Figure |6] shows the cumulative distribution function of 
percentage errors on the weights. This shows that about 80% 
of the stars have less than 20% error on their weights. In all 
further analysis, objects with weight error above 50% were 
discarded. This eliminated 164 objects from the sample, only 
90 of which have {g — r)o > 0.2 and thus belong to the region 
of CMD being analysed. Thus by cleaning just 1% of the 
sample, we remove errant weights without introducing any 
additional bias to our sample . 



3.1.5 Advantages of the adopted weighting method 

Weighting the observed data to allow for parameter- 
dependent incompleteness is a standard requirement in 
statistical sampling, with many types of approaches. The 
method we apply here is both very obvious and very simple, 
and so widely used. To our knowledge it was first applied 
in determina tion of the thick disk abundance distribution in 
section 3.3 of lGilmore et al. 111995), but there may be earlier 

examples. 

In some other investigations (for example, iLarsen et al.l 

(|2010t )). completeness correction is done by comparing the 
dataset in question with a classical power law, which rep- 
resents a smooth fit to the number distribution of stars as 
a function of apparent magnitude, or by comparison to the 
stellar distribution (in colour-magnitude space) predicted by 
a Galaxy model, which might also appropriately take into 
account a model of extinction. The weighting method used 
here also uses the additional parameter of colour. There is a 
definite variation of the value of the weights along the colour 
axis for a given magnitude which makes it important to take 
into account this additional parameter. 

The advantage of this adopted method of accounting 
for completeness is that it is expected to preserve informa- 
tion from any interesting groups that may exist in velocity- 
metallicity space, without being limited by potential mis- 
matches between a model and the actual star counts. 

Similar methods of completeness correction have been 
used by some recent papers analysing stars from the 
Sloan Extension for G alactic Unde r standin g and Explo- 
ration survey (SEGUE). [Cheng et all (|2012bh employ a sim- 
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Figure 5. Observed (dashed blue) and Completeness Corrected 
(solid black) Distribution Functions of [Fe/H] and radial velocity 
for the WGN data fields within 190° <RA< 220°. This covers 
the high latitude region 54° < 6 < 63°, in the inner Galaxy 
312° < I < 347°. The data has been split into fainter stars go > 
18 in the top panels and brighter stars go < 18 in the bottom 
panels. The unweighted metallicity distributions show a slight 
bias to a higher number of metal-poor stars, particularly in the 
faint bins, which are likely to be predominantly halo stars given 
these are high latitude fields. Since these fields are towards the 
inner Galaxy, metal-rich thick disk stars must also be numerous, 
and the weighting method corrects for this, seen more evidently in 
the faint bin. The velocity distributions are in general smoothed 
by the completeness correction. 



ilar method of weighting in g ~ r vs g space to account for 
the incompleteness in old MSTO stars with spectroscopy in 
SEGUE's low galactic latitudes fields. They multiply this 
weight with additional weights - one based on a luminos- 
ity function, to reduce bias against redder metal-rich stars 
in their sample selection; this was found to be mostly con- 
stant, and another based on area covered on the sky to ac- 
count for regions with high extinction that were not consid- 
ered for spectroscopy; th is does not apply to the WGN data. 
iLiu fc van de Veiil (|2012l ) analyse G type SEGUE stars and 
employ weights which are a function of u — g,g — r colours 
and distan ce. Their di s tance s are determined using rela- 
tions from llvezic et al.l l|2008l ) and though they state that 
uncertainties in distance determination will not affect the 
weights, using apparent magnitude instead avoids this po- 
tential problem. 



3.1.6 Comparison to Photometric Metallicities 

llvezic et alj (|2008l ) describe a method of estimating metal- 
licities for those F/G main-sequence stars for which spec- 
tra have not been obtained in SDSS, from their u — g and 
g — r colours, to achieve the same goal as above. Building on 
the classical photom etric metallicity cal i brations applied to 
Galactic structure bv lSandage fc EggenI l| 19591 ), lEggen" et alJ 




20 40 60 80 100 
Percentage Error on Weights 

Figure 6. The cumulative distribution function of the percent- 
age errors on the weights assigned to the stars by our weighting 
method for completeness correction. More than 80% of the stars 
have less than 20% error, confirming that the completeness cor- 
rections derived here are accurate. 



(j 19621 ), IWallersteinI l|l962l ) and lCarnevI l| 19791 ). they first use 
a sample of SDSS stars with spectra to calibrate the metal- 
licities and thereafter use this relation to estimate the metal- 
licities of the complete sample of stars within the same re- 
gion. The photometric metallicity was calibrated using SDSS 
DR6 data in the region Q.2 < g - r < 0.6, 14 < < 19.5. 
The prescr i ption was improved in a follow-up paper by 
iBond et all (|2010l ) (hereafter MWIII) replacing the DR6 
data by DR7. Their calibration is expected to produce good 
results in the metallicity range —2.2 < [Fe/H] < —0.2; at 
the lower limit, [Fe/H] asymptotes due to insensitivity of u- 
g colour to decrease in metallicity, while at the upper end 
there is a systematic underestimation of the [Fe/H] values, 
possibly due to a lack of metal rich stars in the calibration 
data. 

Figure [7] shows a comparison between the metallicity 
distribution functions resulting from WGN metallicities and 
the photometric metallicities of MWIH for the sample field 
RA,Dec=145°,0°. Figure [T^a) shows MWHI estimates for 
those stars which have metallicities directly measured in the 
WGN survey. We can see that at both high and low metal- 
licities the MWIH method shows large deviations from the 
measured values. Figure [7jb) more clearly illustrates that 
MWHI shows systematics of 0.16 dex towards higher metal- 
licities, and very large random errors. The bottom panels 
show how the metallicity distribution functions (MDF) get 
affected as a result. While Figure [TJ^c) shows a compari- 
son between measured WGN metallicities and correspond- 
ing photometric metallicities for the same objects. Figure 
[TJ^d) compares the weighted distribution function in the 
same field with the photometric metallicities for the com- 
plete SDSS photometric sample in that field. Even though 
the distributions produced by the two methods do not differ 
very significantly, the star-by-star comparison shows us that 
the wings of the MWHI distributions are dominated by large 
calibration errors. Photometric metallicities, unsurprisingly, 
blur out the metallicity distribution functions. The photo- 
metric MDF in Figure [Tjd) looks smooth, while the under- 
lying spectroscopic MDF in the same line of sight suggests 
a much more discrete halo plus thick disk two-component 
MDF. 

The weighted dataset that we obtain as a result is thus a 
completeness corrected representation of the WGN selection 
function. 
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Figure 7. Comp arison of WGN m etallicities to photometric 
metallicities from iBond et al. I l l201(t) (referred to as MWIII in 
the plots) for the sample field RA,Dec=145°,0°. [a] Comparing 
the measured WGN [Fe/H] to the [Fe/H] calculated for these 
stars by the photometric metallicity method. The red line shows 
a one-to-one correlation, [b] Difference between the WGN [Fe/H] 
and the photometric [Fe/H]. The photometric metallicities have 
a mean offset of 0.f6 dex towards higher metallicities with re- 
spect to WGN measurements. The la scatter is 0.68 dex. [c] The 
observed metallicity distribution function (MDF) from the WGN 
data compared to the MDF obtained by calculating photometric 
metallicities for the same stars, [d] Completeness corrected WGN 
MDF compared to the MDF resulting from calculating the pho- 
tometric metallicities for the complete SDSS photometric sample 
in that field. 



Section 3.2.3. and verify the robustness of these estimates 
by including the additional constraint of infrared data and 
comparing to estimates from photometric parallax. 



3.2.1 Method 

Distances are estimated by a fitting of the data to the 
main sequence of theoretical stellar isochrones. We start 
with a grid of isochrones with metallicities in the range -2.5 
to 0.0 dex such that the choice of the metaUicity grid evenly 
samples colour-magnitude space. The ages of the isochrones 
are restricted to greater than 10 GYr, since our colour-cut 
selection ensures that only stars from the old Galactic stel- 
lar populations are preferentially sele cted (the thick disk is 
estimated to have an age ~ 12 GYrs iGilmore et al.lll995l '). 
while the halo is similarly old with a most recent e stimate 
of in ner halo field stars to be ~ 11.4 ± 0.7 GYrs l|Kalirail 

Eoii)). 

For each set of photometric and metallicity data de- 
scribing a star, the theoretical isochrones of closest metal- 
licity are first chosen. We then consider those model stars 
on the chosen metallicity isochrone which are within an no- 
colour ellipsoid of the data point. The ncr bracket was set to 
optimise the number of stars for which distances are deter- 
mined. 

Let data point j have k model matches within an no- 
colour ellipsoid. These k matches have the same metallicity 
(by choice of isochrone) and have varying ages greater than 
10 GYrs. The value for a data point j with respect to 
any such model point k is then calculated as follows : 



model^ 



(3) 



where the index x — {u — g,g — r,r — i,i — 2}, each model 
match k having the associated probability given by, 



3.2 Distance Estimates 

We now have a completeness corrected dataset of thick disk- 
halo stars containing information about each star's position 
on the sky and radial velocity apart from ugriz photometry 
and metallicity. In order to use this information to analyse 
the full phase-space structure of the Galaxy we need an esti- 
mate of the third position coordinate - distance. Good stel- 
lar distance estimates are key to making correct inferences 
about global Galactic structure. 

For Galactic distance scales larger than a few hundred 
parsec, while we await Gaia, methods of spectroscopic paral- 
lax are commonly used to estimate distances. Information on 
the temperature and metallicity of the stars is used along 
with stellar evolution models to make an inference about 
the star's absolute magnitude and hence deduce its distance 
given the measured brightness. 

In this work we use a x? technique to fit the data to 
the main sequence of theoretical stellar isochrones in order 
to estimate distances to the stars. In the following sections 
we first outline this method, thereafter describing the use of 
star cluster data to choose the best stellar evolution mod- 
els, and to test the performance of the method. The errors 
on the estimated distances are analysed and calibrated. We 
finally present the distance estimation for the WGN data in 



where m is the number of degrees of freedom. The luminos- 
ity functions appropriate to the isochrones are then used to 
evaluate the probability P(Mfc|agej,) of a star having an ab- 
solute magnitude M^, as suggested by each model match k, 
given its age. 

The probability assigned to each model star k of being 
the best match to the data point j is then given by 

The closest model match to the star j is then found by 
maximising P^^^^i^ , and from this absolute magnitude the 
distance to the star is inferred. 



3.2.2 Testing the distance determination method 

We use star cluster data with well-known distance esti- 
mates to select the set of stellar evolution models which are 
most consistent with the SDSS photometric data on those 
clusters, and thereafter to analyse the performance of the 
method. 
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3.2.2.1 Star Cluster Data lAn et all (|2008l , hereafter 
An08) present photometry for stars in 17 globular clusters 
and 3 open clusters in a value-added catalog of SDSS, de- 
rived by reducing SDSS imaging data in crowded cluster 
fields using the DAOPHOT/ALLFRAME suite of programs. 
Distances were derived to 6 of these globular clu sters and 
2 of t he open clusters in a subsequent paper by lAn et al] 
(2009', hereafter An09) using an isochrone fitting technique 
in the ugriz filter. This reproduced previously available esti- 
mates calibrated from Hipparcos parallaxes in the Johnson- 
Cousins photometric system. 

We applied our method to individual stars in five glob- 
ular clusters looked at by An09 - M5, M3, M13, M92 and 
M15 ranging in metallicity from -1.26 to -2.42 dex and all 
older than 12 GYrs - and compared our resulting distance 
estimates. All these globular clusters have photometry with 
dependable zer o point error f r om A n08. They have [Fe/H] 
estimates from iKraft fc Ivans! (|2003l ) based on high resolu- 
tion spectra. We also chose two open clusters from An08 to 
cover the metal-rich end of the WGN [Fe/H] range - M67 
and NGC2420 with [Fe/H] = 0.00 and -0.37 dex respectively. 
The data for NGC2420 was presented in An08 but not anal- 
ysed by An09. Its distance, metallicity and reddenin g value s 
were adopted by An08 from lAnthonv-Twarog et al] l|2006l l. 

Cluster members were chosen based on a good image 
quality index, a measured characteristic radius close to that 
of the chosen PSF in pixels, and by removing pixels with 
either no detection or saturated detection. The results are 
discussed below, after we consider the optimal choice of 
isochrones. 

3.2.2.2 Selection of isochrones Several sets of 
isochrones are available matched to old stars with SDSS 
photometry. An09, whose distance scale we compare to our 
own, apply the Y ale Rotating Evolutionary Code (YREC; 
ISills et al.l (|2000l ')') isochrones to estimate the star cluster 
distances. They discount the u filter as the models overes- 
timate the flux in the u bandpass causing a 5% difference 
between data and model for colours involving this bandpass. 
More encouragingly, they find that the YREC model colours 
in g — r,g ~ i and g ~ z agree well with the main sequence 
of their globular cluster data. They also adopt an alpha- 
enhancement scheme motivated by the observed amounts of 
alpha elements in field and cluster stars. 

In Figure [5] we comp a re th ree sets of isochrones 
- Dartmouth jpotter et al.l l2008h . YREC and Padova 
ICirardi et al.l I2OO4I ') - with data for the Globular Cluster 
M3. M3 has [Fe/H] = -1.50, estimated age = 13.3 ± 1.4 
GYrs and estimated distance 9.44 ± 1.03 kpc (An09). AU 
isochrones shown have [Fe/H] = -1.5 and age =12.5 GYrs 
(the closest age to that estimated for M3 for which all three 
models have calculated isochrones). An09 use the YREC 
isochrone with [a/Fe] = -1-0.3. The Dartmouth database pro- 
vides isochrones with a coarse grid in alpha-enhancement, 
from which we pick two isochrones of [a/Fe] = and 
-(-0.4. The Padova database does not provide isochrones with 
alpha-enhancement, so we pick [a/Fe] = 0. 

The leftmost panel of Figure [8] plots log luminosity vs 
effective temperature, and this shows that all three mod- 
els are in good agreement with each other, with the excep- 
tion of the systematic deviation of the Padova isochrones 
on the lower main-sequence, which is not relevant here. The 



three isochrone groups employ different methods to convert 
luminosities and effective temperatures to magnitudes and 
colours, which is reflected in the next three panels of Figure 
[S] Here the isochrones have been over-plotted with the data 
for M3, assuming the distance to M3 is as determined in 
An09 (i.e. ~ 9.77 kpc). 

In the Dartmouth Stellar Evolution Model, the model 
colours suffer from inaccuracy in synthetic fluxes for central 
wavelengths bluer than ~ 5000 A thus making the u filter 
very uncertain. Additionally, problems exist in the underly- 
ing SDSS calibration of the u and z filters, which have larger 
statistical errors in their zero points and poorer photomet- 
ric uniformity (both by approximately 50%) as compared to 
the gri filterslj From Figure [8] we see that while the Dart- 
mouth isochrones seem to agree with data in the g — r and 
g — i colours, there is a large deviation (almost 0.5 magni- 
tude near the turn-off of the cluster) m g — z. This suggests 
that there are larger errors associated with the model z fil- 
ter, as compared to gri. The Padova isochrones seem to have 
poorer agreement with data both m g — i and g — z colours. 

In addition to good agreement to data in the gri fil- 
ters, the Dartmouth database provide a finer grid of both 
metallicities (user determined) and ages (0.5 GYrs). Both 
the Padova and YREC databases have more sparse sampling 
at higher ages relevant to the WGN data. Thus we adopt the 
Dartmouth isochrones, excluding the u and z bandpasses 
from the isochrone fitting i.e. in Equation 3 we use only 
colours X = {g — r,r — i}. 

3.2.2.3 Performance The distance determination 
method was tested on the chosen star clusters by picking 
the isochrone of the closest [Fe/H] to the cluster. For all the 
globular clusters, the age of the isochrones was restricted to 
above 12 GYrs. Matches within a 2(j ellipse of two colour 
filters g — r and r — i was used in the fitting. A 2cr el- 
lipsoid was sufficient to obtain distances for at least 80% of 
the cluster data. 

The performance of the distance determination method 
can be analysed by means of an example case of M13 shown 
in Figure [9] The leftmost panel shows a CMD of the clus- 
ter containing all cluster members for which distances were 
determined. The second panel from the left shows a plot of 
g — r vs the absolute g magnitude assigned to the star by 
the distance determination method. The third panel from 
left shows the deviation of this distance estimate from the 
"true" (An09) distance of the cluster as a function of ap- 
parent magnitude g. This panel has been colour coded to 
highlight groups of stars which show extremely deviant dis- 
tances. These deviant groups are as follows: 

• The red points mark the giant stars (i.e stars having 
apparent magnitude g brighter than turn-off) but are incor- 
rectly assigned to be main sequence stars by the routine. As 
a result, their distances are underestimated. These are iso- 
lated using a rough cut of <? < 19 and absolute g magnitude 
>~ 4. 

• The light-blue plus signs on the contrary are main se- 
quence stars (i.e stars having apparent magnitude g fainter 
than turn-off) that have been incorrectly assigned to be gi- 
ants, thus overestimating their distances. These are isolated 

^ http : / / www.sdss .org/dr7 / algorithms / fluxcal . html 
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Figure 8. Wc compare stellar evolution models to star cluster data. Three sets of isochrones - Dartmouth, YREC and Padova - are 
compared with the SDSS photometric data for the globular cluster M3. The adopted properties of M3 are: [Fe/H] = -1.50, age = 13.3itl.4 
GYrs, distance (An ct al. 2009) 9.44 it 1.03 kpc. All isochrones shown have [Fe/H] = -1.5 and age =12.5 GYrs (the closest age to that of 
M3 which has isochrones in all three models). The YREC isochrone has [a/Fe] = -1-0.3; two Dartmouth isochrones with [a/Fe] = and 
-1-0.4 are shown; the Padova isochrone has [a/Fe] = 0. The leftmost panel plots LogL vs Teff, and this shows that all three models are in 
good agreement with each other in the luminosity region of relevance here. The different approaches to convert luminosities and effective 
temperatures to magnitudes and SDSS colours are reflected in the next three panels, where the isochrones have been over-plotted with 
the SDSS photometric data for M3. Prom this comparison we adopt the Dartmouth isochrones restricting their application to g — r and 
g — i colours. 




0.0 0.2 0.4 0.6 0.8 1.0 

g-r 



0.0 0.2 0.4 0.6 0.1 

g-r 



Q 




14 16 18 20 22 24 

g 




5 10 15 
Distance (kpc) 



Figure 9. Results of our distance determination method applied to the globular cluster M13. The left most panel shows the CMD of 
the cluster, the second panel from left shows the absolute g magnitude recovered for each star from the isochrone fitting, the third panel 
from left shows the deviation of estimated distance from the An09 value as a function of apparent magnitude g (red points mark giants 
mis-classified as main sequence stars; blue plus signs mark the main sequence stars which are mis-classified as giants) and the rightmost 
panel shows the distance distribution recovered by the method for all stars (grey dashed line) and for correctly classified main sequence 
stars (solid black line) 



using a rough cut of absolute g magnitude up to the value 
of 4. 

Apart from these anomalies, good estimates of distances 
are assigned to the bulk of the main sequence stars shown 
in black. This can be seen more clearly from the rightmost 
panel in Figure [9] (solid black line) which shows the dis- 
tance distribution recovered for the cluster after eliminating 
all the anomalous star groups (blue and red points) as de- 
scribed above. The blue dashed vertical line which indicates 
the mean distance estimated for the cluster by our method 
matches the green line which indicates the value estimated 
by An09 within 5%. The la spread of the distribution is 



1.1 kpc (placing ~70% of the stars within 15% of the mean 
distance). 

In Figure llOf a) we show the results of testing the dis- 
tance determination method at the metal rich end. All pan- 
els are the same as in Figure |5] There are no giants in the 
cluster, and hence no red points. The blue plus signs mark 
out the main sequence stars mistakenly assigned to be gi- 
ants, as before. But the noticeable difference now is that all 
main sequence stars up to a gr — r colour of 0.6 (which is 
the relevant g — r colour range for the WGN data) seem to 
be mis-assigned as giants. The reason for this can be un- 
derstood from the form of the Luminosity Function for high 
metallicities - for the higher metallicities, the probability of 
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Figure 10. Results of our distance determination method applied to the open cluster NGC 2420 ([Fe/H] = -0.37dex, Age ~ 3.5 GYrs). 
(a) Top panels : Results on applying isochrone fit, restricting ages to above 3.5 GYrs. In this case almost all stars below g — r = 0.6 
are incorrectly assigned to the sub-giant branch, (b) Bottom Panels : Results on applying isochrone fit, restricting ages to above 3.5 
GYrs and additionally restricting the fit to the main sequence region only. Here the problem of mis-classification of main sequence stars 
is resolved. In both (a) and (b), the left most panel shows the CMD of the cluster, the second panel from left shows the absolute g 
magnitude recovered for each star from the isochrone fitting, the third panel from left shows the deviation of estimated distance from 
the An09 value as a function of apparent magnitude g (blue plus signs mark the main sequence stars which are mis-classified as giants) 
and the rightmost panel shows the distance distribution recovered by the method for correctly classified main sequence stars. 



finding a star near turn-off or on the sub-giant branch is 
equal to (or more than) that of finding a star on the Main 
Sequence, as against the lower metallicities where the prob- 
abihty steadily decreases from lower main-sequence going 
towards turn-off. Since we use the luminosity function as an 
additional constraint to determine the absolute magnitude 
assigned to a star, all main sequence stars are identified as 
sub-giants at high metallicities. 

Given the WGN data is comprised mainly of main se- 
quence (mostly turn-off) stars, we can remedy this problem 
by forcing the program to choose model stars that only lie 
on the main sequence (placing a cut-off at turn-off). The 
result of this is illustrated in Figure [TOjb). We immediately 
see that all the main-sequence stars are now assigned good 
distances and the recovered mean distance is within 10% of 
the An09 value. 

While Figures [9l and [TOl illustrate convincingly that the 
method always recovers the correct distance for the cluster 
as a whole, Figure fTHT b) also reveals that the mean recovered 
distance could be a weak function of the photometry itself. 
This trend is further discussed and calibrated in Section 
3.2.2.4. 

Our estimated cluster distances, together with those 
adopted for comparison, are summarised in Table [2] along 
with the various cluster parameters. We see that our method 



reproduces star cluster distances within 10% of the distance 
scale adopted by An09. 

3.2.2.4 Result and Calibration of Errors As noted 
in section 3.2.2.3, while the mean distance recovered for the 
cluster by the distance determination method matches well 
with the true distance, there is a weak trend of the recovered 
distance with photometry. There is also an error which is 
introduced due to the finite [Fe/H] grid that is chosen. In this 
section we analyse and calibrate these different contributions 
to uncertainty in our distance estimates to individual stars. 

(i) Distance error due to photometric errors 
Figure [TT] (top and middle panels) shows the dependence 
of error in distance on error in g — r colour, derived from 
analysis of the photometry for the five globular clusters. We 
calibrate errors based only on the old globular clusters as 
these are more comparable to WGN data in age than are 
the open clusters. Although the photometric uncertainties 
for our field stars are typically less than 3%, for the crowded 
fields of star clusters errors can be larger, allowing sufficient 
dynamic range for this test. On average, 80% of the cluster 
stars have Ug-r less than 0.1. 

The systematic trend in the dispersion of distance er- 
rors at small Ug-r is due to two reasons. The dominant 
effect is the systematic mis-identification of sub-giant stars 
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Table 2. Star Clusters used to test the distance determination method 



Cluster 


[Fe/H] 


Distance 


E{B-V) 


Age 


Dartmouth Isochrone 


Result 


% Error 






( kpc) 




(GYrs) 


[Fc/H] 


Mean Distance 




Globular Clusters 


M5 


-1.26 


7.11 ±0.20 


0.03 


12.2 ± 1.3 


-1.20 


7.36 


4 


M3 


-1.50 


9.77 ±0.32 


0.01 


13.3 ± 1.4 


-1.50 


9.58 


2 


M13 


-1.60 


7.37 ±0.18 


0.02 


14.3 ± 1.1 


-1.60 


7.30 


1 


M92 


-2.38 


8.55 ±0.15 


0.02 


14.4 ±0.9 


-2.40 


8.64 


1 


M15 


-2.42 


10.47 ±0.55 


0.10 


13.5 ± 2.5 


-2.40 


9.87 


5 


Open Clusters 


M67 


0.00 


0.82 


0.03 


3.5 


0.00 


0.84 


2 


NGC2420 


-0.37 


2.51 


0.04 


3.5 


-0.40 


2.23 
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near the turn-off as main-sequence stars below the turn-off, 
which produces systematically offset distances. The second, 
and smaller, of the two effects is the systematic sensitiv- 
ity to age near the turn-off, with the method showing some 
bias to younger isochrones weighted by a higher probabil- 
ity from the luminosity function. Not surprisingly, restrict- 
ing the isochrone search to the age adopted for the parent 
cluster diminished the distance bias. Of course, a single age 
cannot be applied to the mixed population field star WGN 
data, for which only a reasonable lower limit of age can be 
imposed. This distance bias, though, defines our error limits, 
and is sufficiently small as not to limit our later analysis. 

(ii) Distance error due to uncertainty in [Fe/H] 

In our method, we have assumed that the [Fe/H] estimate 
is accurate, and hence choose the isochrone with metallicity 
closest in value to this measured value. Here we quantify 
the uncertainty in the distance introduced by uncertainties 
in metallicity , which we determined to be ~ 0.2 dex in 
[Fe/H]. We do this by applying the distance determination 
method to the stars in the test clusters while varying the 
[Fe/H] of the chosen isochrone around the best value. The 
bottom panel of Figure [11] shows the results in the case of 
four of the globular clusters (M15 is not shown here as it has 
same [Fe/H] as M92). The WGN and SDSS data,to which 
our method will be applied, have [Fe/H] accuracies ~ 0.2 
dex, and from the figure we determine that the systematic 
error on distance due to 0.2 dex uncertainty in [Fe/H] is 
within 10%. 

(iii) Distance error due to luminosity mis-classification of 
stars 

By restricting the isochrone fit only to the main sequence 
region, we ensure that no main sequence star is mis-classified 
as a giant (as illustrated in Figure fTOT bl'). But the sub-giants 
in the WGN data near turn-off are mostly (more than 95%) 
mis-classified as main sequence stars (as illustrated by the 
red points in Figure [9|. From point (i) in this section we 
have shown that the there is a systematic error due to this 
mis-classification of sub-giants near the turn-off of within 
10%. From Figure [9] we note that for sub-giants redder than 
turn off the distance errors can blow up to as large as 50%. 

We also note that in the case of clusters M5 and M92, 
there is no upward trend in errors near small (Jg-r, as com- 
pared to the noticeable trend in M3,M13 and M15. This is 



because of the presence of fewer giants (less that 0.5%) in 
the data of the former two globular clusters, as compared to 
the latter in which more than 4.5% of the cluster members 
are giants. 

The fraction of a field star sample which will be sub- 
giants rather than MSTO stars was a popular topic of dis- 
cussion in the 1980's star-count lite rature. Some of man y 
papers discuss i ng th e issue in c lude [Bahcall et al.l (1 19851) . 
'Gil more etahl l| 19951 '). iFenkartI (ligsgV and iGilmore 



(198^ All are consistent with results that can now be read- 
ily reproduced using the on-line Besangon star count model 
as shown in Figure [T2l This confirms that sub-giant contam- 
ination in a colour and magnitude selected survey of thick 
disk and halo stars with photometric criteria like those of 
this WGN sample will be typically at the 5% to 10% level, 
increasing to 13-15% in some lines of sight. 

Thus we can say that sub-giants, certainly present in the 
WGN data, will be mis-classified as main sequence turnoff 
stars, introducing a systematic underestimation of their dis- 
tance by 20 to 50%. We expect this contamination from 
sub-giants to be small (below 10% per field) and thus not 
affect our present analysis in a significant way. 

3.2.3 Distance estimates for the WGN data 

We thus derive distances for stars in the WGN dataset via 
the isochrone fitting method described above using g — r 
and r — i colours. We use a set of 17 i sochrones from the 
Dartmouth Stellar Evolution Database (|Dotter et al.ll2008l 'l 
for metallicities in the range -2.5 to 0.0 dex (and no alpha- 
enhancement), such that the choice of isochrones evenly 
samples colour-magnitude space. We fit to only the main 
sequence of the isochrones and restrict their ages to greater 
than lOGYrs. The corresponding luminosity functions (age 
and metallicity dependent) were adopted from the same 
database. A 4(t colour bracket was required to find matches 
for lower main sequence stars with small (Jg~r where the 
isochrone sampling is sparse. Distances could be estimated 
for 8965 objects out of the 9404 in the sample. The few 
hundred stars for which the routine was unable to deter- 
mine distances are those which do not fall on the thick 
disk/halo CMDs either because they are extreme blue stars 
((? — r < 0.2) or because they are likely binaries, variables or 
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Figure 11. Errors in distance estimation, as a function of un- 
certainties in photometry and metallicity, are plotted for five star 
clusters (M5, M3, M13, M92 and M15). Top panel; the disper- 
sion, and Middle panel: the percentage error in distance due to 
errors in the single-star photometry. The bias in the dispersion 
at small cTg—r is due to the mis-identification of sub-giant stars 
near the turn-off as main-sequence stars. Lower panel: The error 
in distance due to a systematic error in metallicity, determined by 
applying the distance determination method to the star cluster 
photometry while adopting an isochrone with [Fe/H] offset from 
the notional cluster value. The systematic distance error corre- 
sponding to 0.2dex uncertainty for [Fe/H] is less than 10%. 



have mis-calibrated metallicities. The prior assumed for the 
age (older than lOGYrs) also produces a bias against metal- 
rich stars. The distance distribution peaks at just above 2 
kpc, most stars being within a heliocentric distance of 6 kpc. 

We now consider adding additional photometric infor- 
mation, as a further test of distance scale systematics. 

3.2.3.1 Including infrared photometry - testing 
systematic distance bias In order to further check the 
robustness of our distance estimates, particularly against 
systematic bias, we tested the value of incorporating infrared 
photometry into the isochrone fitting, combining J-, H- and 
K-band data with the gri-h&nd data. 

UKIDSS may be considered as th e near-infrared coun - 
terpart of the SDSS Imaging survey l|Hewett et al.ll200d ). 
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Figure 12. Percentage of sub-giants (Grey Diamonds) and thin 
disk (Black asterisk) stars in each field as produced by the Be- 
sangon Model of the Galaxy on application of the WGN selection 
function. This gives an estimate of the expected contamination 
from subgiants and thin disk for a colour and magnitude selected 
survey of thick disk and halo stars with photometric criteria like 
those of this WGN sample. The plot confirms that sub-giant con- 
tamination in the WGN sample will be typically at the 5% to 
10% level, increasing to 13-15% in some lines of sight (discussed 
in Section 3.2.2.4). It also shows that thin disk contamination is 
expected to be minimal in the WGN sample (discussed in Sec- 
tion 4.3), ranging between 5 and 15% in the low latitude fields 
(i.e near Ra~ 133°,Ra~ 249° and Ra~ 313°). 



The Large Area Survey (LAS) covers 4000 sq. deg. at high 
Galactic latitudes, with imaging in YJHK filters to a depth 
of K=18.4. This coverage includes a LAS equatorial block 
which overlaps with most of the sky coverage of the WGN 
survey. 

We cross-matched the WGN dataset with the UKIDSS 
LAS data. This provided 3409 matches out of 9404 stars in 
WGN. We used the JHK filters in tandem with the gri filters 
to estimate the distances by our isochrone fitting method i.e. 
we now use x = g ~ r,r ~ i, J — H, J — K m Equation 3. The 
distance determination proceeded as described above. We 
compare the new distance estimates, obtained using JHK 
both independently and with gri, referred to as Djhk and 
DgrijHK, to the distance estimates from just gri (Dgri) in 
Figure [T3l 

We see that Djhk shows a mean offset of 0.08 kpc with 
respect to Dgri, though with a large scatter of la — 0.48 
kpc. This greatly improves in the case of Dgrijun where 
the mean offset reduces to -0.03 kpc and a much tighter 
spread of la = 0.17 kpc This shows that there is no bias in 
the estimated distances using gri relative to the combined 
optical-infrared approach. In this same figure, we also illus- 
trate the effect of retaining the z filter in the left most panel. 
Dgriz has a systematic bias, being offset by around 0.2 kpc 
from Dgri- 

Thus by comparing our distance estimates with those 
obtained using tighter constraints (i.e. using JHK along 
with gri) , we confirm that there is no systematic bias in our 
distances. 

3.2.3.2 Comparing isochrone-fit distances with 
photometric parallax distances We also compared our 
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Figure 13. This figure compares isochrone-fit distance estimates using optical only, infrared only, and combined optical 
plus infrared photometry. The first three panels plots the relative error in distance estimates from using griz filters (blue), 
the UKIDSS infrared JHK filters (green) and by combining the gri and JHK filters (red) with respect to distance estimates 
using gri filters {Dgri). Dgriz has a systematic bias, being offset by around 0.2 kpc from Dg^i- Djuk shows a mean offset 
of 0.08 kpc with respect to Dgri, though with a large scatter of la = 0.48 kpc. This greatly improves in the case of DgrijHK 
where the mean offset reduces to -0.03 kpc and a much tighter spread of Icr = 0.17 kpc. 
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Figure 14. We compare the photometric parallax method of estimating distances to the isochrone fitting technique used 
in this work. In the first three panels distance estimates fro m our isochr o ne fit ting me thod are compared to distances 
from photometric parallaxes, us ing the recipes prescribed in llvezic et al.l ||2008| ) (red), ISmith et al.l ||2009| ) (blue), and 
ISmith, Whiteoak fc Evanj {Soli) (orange) in terms of relative distance error. 108 systematically overestimates the distances 
with an average offset of ~ 0.09 kpc, showing larger deviations at larger distances. The match improves in the case of S09 
with a mean offset of ~ 0.08 kpc, while the SWE12 method matches the distances computed in this work most closely. The 
mean offset with respect to SWE12 is ~ —0.04 kpc, with a spread of ~ 0.12 kpc. 



distances to those obtain ed by applyins the method of pho- 
tometric paralla:x used bv llvezic et al.l l|2008l 'l (hereafter 108). 
In this method, an empirical relationship between absolute 
magnitude and colour is first adopted, which is indepen- 
dent of both age and metallicity and only describes the 
shape of the colour-magnitude sequence below turn-off. This 
zeroth-order relationship is usually adopted by finding the 
best-fit function to the CMD of a globular cluster of known 
abundance and whose distance is known fairly accurately by 
means of some other reliable method of distance measure- 
ment. Then correction terms are introduced which shift this 
relation to represent the true age and metallicity. 

108 perform an empirical fit simultaneously to the pho- 
tometry of stars of five globular clusters from SDSS to obtain 
the basic photometric parallax relation. They then use the 
M13 globular cluster to determine the correction term due 
to the effect of age on the CMD. This age correction is cal- 



culated under the constraint of a fixed value of metallicity - 
which in their case is the mean metallicity of M13 ([Fe/H] = 
-1.54), which is close to the mean metallicity of the Galactic 
halo. Strictly speaking, this will not hold for thick disk stars 
whose mean metallicity, at least in local samples, differs by 
the order of +1.0 dex from that of halo stars, even if the 
two populations have approximately the same age, within 
an uncertainty of 2 GYr. 108 adopt this correction on the 
basis that this will produce an error of less than 0.2 mag 
for disk-like stars (which corresp onds to an error in d istance 
estimate of 1 kpc in our case). ISmith et al.l (|2009| ) (here- 
after S09) propose an improved correction which takes into 
account that the morphology of the colour-magnitude re- 
lation near turn-off depe nds significantly on both age and 
metallicity. More recentlv fSmith. Whiteoak fc Evanj (|2012l ') 
(hereafter SWE12) modified this to construct a global turn- 
off correction applicable to a mixture of disk and halo stars. 
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All the above methods were applicable in the colour region 
of 0.3 < g ^ i < 0.6. This colour restriction eliminates close 
to 2000 stars from the WGN dataset, the majority of these 
being thick disk turn-off stars with g — r ^ 0.6. 

Selecting only those stars from the WGN dataset which 
fall in the photometric calibration g — i colour range, the 
distances computed here using our isochrone matching tech- 
nique may be compared to the values obtained by apply- 
ing the methods of 108, S09 and SWE12. The comparison 
is shown in Figure [ill 108 systematically overestimates the 
distances to the stars, with larger deviations at larger dis- 
tances. The 108 values are on average offset from distances 
from this work by ~ 0.09 kpc. The match improves in case 
of S09 with a mean offset of ~ 0.08 kpc, while the SWE12 
method match the distances computed in this work most 
closely. The mean offset with respect to SWE12 is ~ —0.04 
kpc, with a spread of ~ 0.12 kpc. 

This provides a final consistency check on our distance 
estimates. 

The method of photometric parallax works well in the 
case of faint main-sequence stars and single stellar popula- 
tion, but is limited by the availability of calibration stars 
over a large range of stellar mass and metallicity. In the 
case of a mixture of populations, the age correction in the 
photometric parallax method requires an assumption about 
the fraction of stars from different populations which S09, 
for example, determine based on a kinematic and chemical 
segregation of stars. Our method based on stellar evolution 
models avoids any such assumption and is more robust in 
accounting for the age and metallicity effects for turn-off 
stars. 

We thus have estimated distances to 8965 stars in the 
WGN dataset, with 90% of them having uncertainty within 
10%. 

Figure [15] shows the spatial coverage of the data. The 
data covers a range of Galactocentric radii (projected dis- 
tance on the plane) between 4 and 12 kpc. There is also 
coverage above and below the Galactic plane, extending up 
to a height of 6 kpc. This provides the scope to investigate 
radial and vertical gradients in the abundance and kinematic 
parameters, as well as to investigate spatial distributions in 
a coherent way. The coverage of the data on the Galactic 
plane, shown by the X vs Y plot, shows the potential of the 
data to investigate symmetries in the Galaxy - about the 
plane, about the cardinal direction 1~ 270° and about the 
centre of the Galaxy. 

3.3 Proper Motions and Global Space Velocities 

Proper motions (along RA and Dec) for the WGN stars 
are available from SDSS DR7. Using these along with the 
distance estimates obtained as described in the previous sec- 
tion and the measured radial velocities, we can obtain right- 
handed Galactic U,V,W velocit i es usi ng the conversion out- 
lined in [Johnson fc SoderblomI l|l987f ). These are then con- 
verted to cylindrical coordinates Vr, V^, Vz by simply rotat- 
ing the coordinate system about the W axis. 

In the cylindrical coordinate system, the contribution 
of the three velocity components to the line of sight (LOS) 
velocity at any point in space specified by the Galactic lon- 
gitude I, Galactic latitude b and distance from the sun D are 
given as 



Vlos = ArVr + A^V4. + AzVz (5) 
where 

Ar — sin(Q)cos(b) (6) 
A,/, — cos(Q)cos(fo) 
Az ~ sin(fo) 

and where a is the angle between and the line of sight 
direction which depends on l,b and D and is given as follows: 



Q = TT — / — 6gc 
where 6gc 



(7) 



tan 



73cos(b)sin(l) 



Ro — I?cos(b)cos(l) 

and Ro is assumed to be 8.5 kpc l|Kerr fc Lvnden-Bell[[l98^ ). 

Here, Vlos is the line of sight velocity component in the 
rest frame of the Galaxy which is related to the measured 
heliocentric radial velocity Vhet by 



Vlos ~ Vhel + Vcorr + VlSRcot 



(8) 



where Vcorr corrects for the peculiar motion of the 
sun [we adopt (£/(:) Va,W(T)) = (11.1,12.2,7.2) km/s 
(Sch onrich et al.[[2010l )] with respect to the Local Standard 
of Rest 

Vcorr ~ (70cos(/)cos(fe) + V0sin(/)cos (6) -I- W(Dsm{b) (9) 

and VLSRcorr converts the measured quantity into the rest 
frame of the Galaxy: 

VLSRcorr = VLSi{sin(Z)cos(6) (10) 

Vlsr has been taken to be 220 km/s ([Kerr fc Lvnden-Bell[ 
[ 19861 ) . 

Given the spatial distribution of stars in our dataset, 
the relative contributions of the three Galactic Space ve- 
locity components to the measured line of sight velocity 
varies widely across the fields. The contribution of V^ goes 
to zero in the direction towards the Galactic Centre, while 
the Vr contribution goes to zero perpendicular to this direc- 
tion {I ~ 90° and 270°). This is reflected in the variation of 
uncertainties in the velocities of each of these components 
along Galactic longitude. While Vr is most accurately de- 
termined towards the Galactic centre from Vlos alone, V^ 
is more reliably determined in the perpendicular direction. 

In this work, we are particularly interested in the global 
behaviour of the V^ parameter, as this is a primary thick 
disk-halo discriminator. Since the sensitivity of the line of 
sight velocity to varies considerably across our flelds, 
large errors in proper motion can significantly distort the 
true V^ distribution. Figure [16] analyses the contribution of 
errors from the two different sources, namely from distance 
and proper motion errors, to the V^ distribution. Considered 
in this figure are only fields centred on the cardinal direc- 
tion / — 270° , which makes the observed radial velocity 
more sensitive to V^ as compared to Vr. The leftmost figure 
shows the distribution of the computed if both distance 
and proper motion measurements are considered to be ac- 
curate. The middle and rightmost plots show how the Prob- 
ability Distribution Function (PDF) is altered if there is a 
contribution of error from only distances and only proper 
motions respectively. The typical errors in proper motion 
are around 3mas/yr while typical errors in distance are 0.2 
kpc. When considering errors only from proper motions (i.e. 
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Figure 15. The spatial extent of the WGN survey data in XYZR galactic coordinates. In both plots, the black cross shows the position 
of the Sun. Colours indicate the number of stars in each spatial bin, colour coded as shown. 
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Figure 16. Space motion error analysis for fields in the direction I = 270°. Left Panel: distribution derived using radial velocity, 
proper motions and distance, while assuming that there are zero errors on distance and proper motions. Middle panel: V0 distribution 
when convolved with the contribution of error from distances alone (i.e. setting = 0). Distance errors are less than 10%, thus 

ensuring that 80% stars have less than 10% error on V^. Right panel: distribution when convolved with the errors from proper motions 
alone (i.e. = 0) which are typically 3mas/yr resulting in 85% stars having greater than 20% error on V^. Evidently proper motion 
uncertainties dominate the error budget. 



setting an = 0) only 15% of the stars have ay^ within 10%. 
The typical errors on are ~ 13%. Almost 40% of the 
stars have errors on greater than 20% thus smearing out 
the distribution on the rightmost panel of Figure [161 In con- 
treist, when considering errors only from distances (setting 
""(Jo 5 ""Mi ~ 0) 9i0% of the stars have ay^ better than 10% 
with typical errors of 3%. 

This clearly illustrates that proper motions are the chief 
contributors to errors in the velocity components, a factor 
we consider in the analysis below. 



4 ANALYSIS 

We now proceed to analyse the WGN survey data, consider- 
ing in turn the validity and numerical values of single para- 
metric descriptions of the thick disk, for kinematics and spa- 
tial scale length. We then consider, and discover, spatially 
localised substructure in the survey. 



4.1 Kinematic and Number-Density Modelling in 
Observable Space 

When analysing a set of distances, abundances and kine- 
matics for a mix of stellar populations one approach is to 
minimise error propagation by looking for models in phys- 
ical space that most closely reproduce probability distribu- 
tion functions in observable space. In our sample, the radial 
velocities are of higher internal precision than are the other 
contributions to space motions (cf. Figure [T6|l . so we proceed 
accordingly. 

We model the data in the one dimensional space of line 
of sight velocity Vlos- The data are assumed to be com- 
pletely described by a two-component model consisting of 
the thick disk (denoted as TD) and the halo, using the dis- 
cussion associated with the selection function earlier to jus- 
tify exclusion of thin disk populations. Initially, the three 
cylindrical velocity components i.e. Vr,V4,,Vz are assumed 
to be distributed as Gaussians for the thick disk and halo 
independently. 

A maximum likelihood method is employed to deter- 
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mine the best model given the data. The hkelihood of a 
particular model <1> given n data points X — xi, X2..., x„ is 
then given by 

n 

L{$\X)^l[p{x.\$) (11) 

1 

where 

P(x,|$) = /P(a;.|$™) + (1 - 

= fP{^'^°{VLOS,^\VLOs)) 

+ (1 - f)P{^''''"'{VLOS,.\VLOs)) 

with / being the normalisation of the thick disk component 
with respect to the total number of stars. 

Since we do not know the shape of the distribution of 
Vlos we can recast this in terms of Vr, V4,, Vz which have 
(assumed) Gaussian distributions, and then marginalise over 
the unkno wn paramet e rs. T his is done using the method 
outlined in lSirko et al.l (|2004l ) in their Appendix A. 

The probability of observing a star with a line of sight 
velocity between Vlos and (IVlos is given by 

/ <t}R{yR)<t,4,{V4,)4,z{Vz)dVRdV^dVz 

-00 J —00 

with 

^.(^.)rfK = y^exp(^(i;^-^). (12) 

Using the form of Vlos given by equation O equation 
[12] can be analytically solved to obtain the following simple 
solution 

4>{yLos) = ^=eM~^] (13) 

where 

S = OrA^ + cr^A^ + o\A%, 
A = —Vlos + Mfl^fl + ^-4,^4, + ^izAz 

with Ar, A^, Az being defined in Equation 6. 

The best relative normalisation / is obtained indepen- 
dently for each field for any given configuration of all the 
other four free parameters. 

4.2 Maximum likelihood determination of the 
thick disk kinematic parameters 

In this analysis we seek first to find a maximum likelihood 
description of the velocity components of the thick disk. The 
free parameters in the model are fiv^ (henceforth written 
simply as V^) and av^j , fv^ , ctv^ of the thick disk, and /, 
the relative number of thick disk stars to halo stars for each 
field. All kinematic parameters describing the halo, together 
with the mean Z and R velocities of the thick disk, are held 
as fixed parameters in the model. The mean Z and R veloc- 
ities of the thick disk and mean motions of the halo are set 
to zero, assuming equilibrium and that there is no net rota- 
tion in halo. All those ~ 8900 stars with good photometry, 
Cg-r ^0.1 and weight error less than 50% were included 
in the modelling. To ensure kinematic parameters appropri- 
ate to SDSS survey results covering comparable volumes to 
those investigated here, we adopt the values of the velocity 



disper sions of the inner halo determined by ICarollo et al.l 
(l2010l ). These values are summarised in Table [S] 

Figure fTTl shows the Galactocentric velocity distribution 
for the data and the corresponding distribution produced by 
the best fit model resulting from our ma^ximum likelihood 
analysis, where the best fit values were found by scanning 
through parameter space with a step size of fi, a— 1,1 km/s. 
In this plot, the number of stars is normalised to be the 
same in each field in order to pick out features easily by eye. 
The best fit parameters recovered by the Maximum Likeli- 
hood modelling were V^, ay^ — 172, 49 (km/s) and avpt , crvz 
= 51, 40 (km/s). The rotational velocity is only marginally 
lower than a canon ical thick disk rotation (~ 180 km/s, 
ICarollo et all (|2010l )l. The contour plots localising maxima 
and indicating error bounds and parameter correlations for 
the various parameters are shown in Figure 1181 These fig- 
ures show that this simple kinematic model is a good global 
description of the data, while our best-fit kinematic parame- 
ters are in go o d agr e ement with earlier det erminations (e.g. 
ICarollo et ahl (|2010l ). IChiba fc Beerd (|2000l 'l). 

The significance of the difference between model and 
data (i.e difference between model and data divided by the 
Poisson error in number counts in each pixel) is plotted in 
the bottom-most panel of Figure [T71 A significant mismatch, 
which we define as several adjacent pixels each of which is 
more than 2a from the model with the same sign, is evident 
only near RA ~ 313°. Here many stars seem to have higher 
line of sight velocities than predicted by the model. 

The region of mismatch covers 312° ^ RA ^ 318°. It 
spans low, southern latitudes —30° ^ & ^ —26° and lies 
in the inner Galaxy 48° I ^ 50°. It has radial velocities 
~20km/s higher than do the typical thick disk stars in these 
lines of sight. Hereafter, we will refer to this region as R313. 

4.3 Sub-structure in the Galactic thick disc? 

Is the mismatch seen between the global best-fit kinematic 
model and the data due to substructure in the thick disk, or 
is it just an artefact of incorrect extinction or systematic er- 
rors in distance? Are the stars with high radial velocity with 
respect to the thick disk just high contamination from thin 
disk stars at these low latitudes? We analyse the possibilities 
in this section. 

R313, being at the lower Galactic latitude range of the 
WGN fields that looks into the inner parts of the Galaxy, is 
likely to have a larger contamination of thin disk stars. The 
high rotational velocity of the thin disk will show itself with 
radial velocities higher than that of thick disk by similar 
amount as seen in this direction. To investigate whether an 
inadequate model of thin disk stars contributes to the mis- 
match, one may consider comparing R313 with the corre- 
spondingly low Galactic latitude WGN field near RA~ 133° 
{£ — 228°, fe — -1-27°). We found, especially by comparing the 
metallicity distributions in these two directions, that there is 
no significant contribution from a thin disk population. By 
looking at the number of stars in each field as given by the 
Besangon Model of the Galaxy (iRobin et al.ll2003 ) on ap- 
plying the WGN colour-magnitude cuts, we verify that the 
number of thin disk stars that pass the selection function is 
minimal , the contamination being close to in most of the 
high latitude fields and between 5 and 15% in all the lowest 
latitude fields, as shown in Figure 1121 The thin disk kine- 
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Table 3. Kinematic parameters adopted and determined in the Galactic model 



Parameter Value fixed (km/s) Value determined (km/s) Reference 



Halo 

Vr,V4„Vz 

<^Vr 



0, 0, 
150 ± 2 
95 ± 2 

85 ± 1 



assumed 

CaroUo et al. (2010) 
CaroUo et al. (2010) 
CaroUo et al. (2010) 



Thick Disk 
Vr, Vz 

avz 



0, 



172 ± 2, 49 ± 2 
51 ± 3 
40 ± 3 



assumed 
this study 
this study 
this study 
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RA 

Figure 17. Galactocentric velocity vs RA is plotted for the weighted data (top panel) and for the model with our best fit thick disk 
rotational parameters of V^,a-v^ = 172,49km/s and (JVn^'^Vz = 51,40km/s (middle panel). For these two panels, colours indicate the 
number of stars in each bin, as indicated. The dashed lines correspond to a rotational velocity of 172 km/s, the determined best fit mean, 
at 2 kpc. Both data and model plots have been normalised to have 1000 stars in each field in order to bring out features more clearly. All 
the velocity parameters of the halo are fixed in the model, (see Table [Sj . The relative normalisation between thick disk and halo for each 
field is allowed to be a free parameter. Lower panel: this shows the significance of the difference between the model and the data. In this 
panel, blue to green denotes positive deviation ranging from to +3cr while red to orange denotes to —6a. The dashed lines in the lower 
panel correspond to a rotational velocity of 172 km/s and 100 km/s at a distance of 2 kpc, and are illustrative. Galactic latitude is also 
indicated at the top of upper panel. 
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Figure 18. The likelihood contours of the four thick disk kinematic parameters determined in the maximum likelihood fitting of the model 
- V^,a^, an, crz- The red box in each case marks the value for which the likelihood is maximised. 



matic signal, if present, may be unidentifiable towards the 
Galactic centre, near fields RA~ 249°, but RA~ 133° being 
similarly inclined to the centre-anti-centre line as RA~ 313°, 
we would expect the thin disk stars to equally show them- 
selves in the line of sight velocity distributions in these fields 
as well, but this is not the case. The data match our kine- 
matic model well near RA ~ 133°. Hence we conclude that 
a misrepresentation of the high-latitude symmetric thin disk 
cannot be the cause of the excess numbers of stars with high 
line of sight velocities. For similar reasons, the fact that we 
may have ignored any vertical velocity gradient in the kine- 
matic model for the thick disk is an unlikely cause of the 
mismatch. 

Might incorrect interstellar extinction be at fault? If 
one underestimates the extinction one overestimates the dis- 
tance, which can make stars in this region appear to be 
moving faster. The estimated distance of stars in R313 is on 
average between 2 to 3 kpc. However, even a 25% distance 
error increases the Vgsr at RA ~ 313° only by 5 km/s for 
a rotational velocity of 180 km/s. Therefore, overestimated 



distances due to incorrect extinction cannot produce as sig- 
nificant a shift in Vgsr as is being observed here. This would 
also not explain the excess star counts discussed in the next 
section. 

The above checks confirm that the velocity mismatch 
seen in region R313 is real and due to some substructure. 

4.4 Modelling number counts - the thick disk 
scale length 

We look at how many stars in the R313 fields contribute 
to the velocity mismatch. For this, we first assume that the 
substructure is disk-like with similar Vr and Vz distribu- 
tions as the thick disk, but a higher mean V^. We fix the 
halo kinematics as before and fix the thick disk kinematics 
to the values determined from the global fit in Section 4.2. 
We isolate four WGN fields that cover region R313 and fit 
for and cr^ of the substructure. We obtained a maximum 
likelihood fit of V^,a^ =206,36 km/s. As shown in Figure 
[19] the substructure has a 35% and 30% contribution to the 
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Figure 19. The panels show hne of sight velocity for four fields near R313 where the mismatch in velocities between data and model 
is apparent. Data in these fields are simultaneously fit by three components - thick disk, halo and a third component representing the 
overdensity. Halo kinematics are held fixed as before (Table [Sjl. Thick disk kinematics are adopted from the fit result in Section 4.2 i.e. 
Vij,,(TY^,(TYii,'^Vz ~ 172,49,51,40 km/s with mean Vji,Vz fixed to zero. The overdensity is assumed to be thick disk like, albeit with 
a different rotation velocity. All parameters for the overdensity are fixed to thick disk values, except , crv^ which are fit by maximum 
likelihood. The best fit value obtained was 206,36 km/s. The histograms show the unweighted data, the blue curve show the thick disk 
model, red curve shows the overdensity model and the black curve shows the total model in each case. The RA of the field is indicated 
above each panel. The contribution of this extra component is as high as 35% and 30% in the left two panels while it drops to 5% in the 
last two. 



fields centred on RA 313° and 315° respectively while the 
contribution falls sharply to 5% in the fields centred on RA 
317° and 319°. We should expect these stars to show up as 
an excess in number count distributions of the thick disk. 

We therefore go on to analyse the number of stars per 
field present in the SDSS DR7 data adopting the same 
colour-magnitude cuts as the WGN selection function, and 
compare this to the number count predictions given by the 
Besangon Model of the Galaxy. This comparison is shown in 
Figure 1201 where significant differences between the model 
and the data are evident in the RA ranges which correspond 
to the lowest Galactic latitudes (RA~133° and RA~313°). 
The default extinction in the Besangon model is smaller than 
in SDSS data, so we also show model predictions adopting 
the SDSS DR7 mean extinction in each field, confirming a 
real difference between model and data. 

In order to understand the reasons for this difference be- 
tween model and data, we describe here a simple toy-model 
analysis. We generated the star count density distribution of 
the thick disk and halo using respectively a double exponen- 
tial disk parametrised by scale height hz and scale length 
Hr to describe the density distribution of the thick disk and 
a flattened power-law with fla ttening q = 0.6 and index 
Q = 2.3 to describe the halo IjDeason. Belokurov fc Evand 
[soil), with the number of halo to thick disk stars fixed 
within the Galaxy. In Figure [21] we show the number dis- 
tribution generated by this model for the WGN fields and 
selection function. 

In the top panel of Figure [211 fez is kept fixed at 0.9 kpc 
and Hr varied from 2 to 8 kpc, while in the middle panel Hr 
is fixed at 4.0 kpc while hz is varied from 0.7 to 1.3 kpc. From 
the plots, we can see that the number counts at lower Galac- 
tic latitudes are the most sensitive to thick disk structural 
parameters - as the thick disk scale length is increased (or, 
with less effect, the scale height is decreased), the number 
of stars near RA~135° and RA~315° increases. The data 
are mostly sensitive to the thick disk radial scale length and 
a sufficient range of fields has been targeted to provide use- 
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Figure 20. The surface number density distribution of stars per 
WGN field is compared between SDSS DR7 and the Besangon 
Model of the Galaxy for a magnitude range of f6 ^ g 17.5 and 
same colour cuts as the WGN selection function, black : SDSS 
DR7 data ; orange : Besangon Model values ; green : Besangon 
model number counts corrected with the corresponding extinction 
for each field from SDSS DR7. SDSS data shows an excess of 
stars in comparison to the Besangon model in two regions - RA 
less than 140° and RA in the range 313° to 330°, which figure [211 
shows is the region sensitive to thick disk structural parameters. 



ful constraints on those structural parameters. This analysis 
thus illustrates the utility of the range of fields selected for 
observation by WGN as a test of thick disk structural pa- 
rameters. 

We define two figures of merit, Ratiol, being the ratio 
of the number of stars per square degree selected by the 
WGN criteria in the field at RA=313° to the corresponding 
number for the field at RA=249°; and Ratio2, the corre- 
sponding number comparing the star counts in the field at 
RA=133° to those in the field at RA=249°. These figures 
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Figure 21. A simple illustration of the sensitivity of star counts 
in the WGN fields to the structural parameters of the thick disk. 
The figure shows number count distributions produced by a sim- 
ple density model for thick disk-halo : a double exponential disk 
parametrised by scale height hz and scale length hji for the thick 
disk and a flattened power-law with flattening q = 0.6 and index 
a = 2.3 for the halo. Distributions are obtained for the WGN 
lines of sight for stars in a 6 kpc region centred on the Sun. The 
global ratio of halo to thick disk stars in the Galaxy is held fixed. 
Top: hz fixed at 0.9 kpc. hji varied from 2 to 8 kpc in steps of 
2 kpc. Middle: hn fixed at 4 kpc. hz varied from 0.7 to 1.3 kpc 
in steps of 0.2 kpc. We can see that increasing hji or decreas- 
ing hz significantly changes the number counts near RA~135° 
and RA~315°, which correspond to the lowest Galactic latitudes 
probed. Bottom: Comparison of the simple parametric star count 
galaxy model with star count distributions in the WGN fields 
predicted by the Besangon model. The values of the two figures 
of merit defined in the text, Ratiol and Ratio2, are given for the 
distributions for quantitative comparison. The model with hji, hz 
= 2.0,0.8 kpc (blue) reproduces the number distribution of the 
Besangon model (green) adequately, showing these are equivalent 
exponential scale lengths to the model for present comparative 
purposes. 



of merit quantify how the change in the thick disk struc- 
tural parameters hn , hz affect the shape of the star counts 
distribution as a function line of sight. The star count distri- 
butions produced by the Besangon model has Ratiol=0.66 
and Ratio2=0.20 while the observed distribution from SDSS 
data has Ratiol=1.69 and Ratio2=0.37. This quantifies the 
discrepancy between data and model evident in Figure [20l 

By experiment we see that {hu,hz) = 2.1,1.0 kpc re- 
produces the Besangon number distribution well (shown in 
bottom panel of Figure I21|l . That is, the Besangon model 
description of the thick disk, for the observational param- 
eters of relevance here, may be described as a double ex- 
ponential with (hit,hz ) = 2.1,1.0 k pc. T his is consistent 
with values quoted in iRobin et al.l (|2003l ) as used in the 
model for the thick disk density profile which is a modified 
exponential having scale length of 2.5±0.5 kpc, and scale 
height of 0.8±0.5 kpc. As shown in Figure [20l however, that 
model substantially and systematically under-predicts the 
observed SDSS star counts. Increasing the radial thick disk 
scale length to the value ha = 4.1 kpc, as shown in Fig- 
ure [221 makes the shape of the number distribution consis- 
tent with SDSS d ata.This number i s in agreement with the 
value obtained bv ljuric et al.l l|2008h who analyse M dwarfs 
in the solar neighbourhood (at distances less than 2 kpc) 
from SDSS data to estimate a scale length of 3.6 kpc for 
the thick disk. While most previous works which estimate 
the scale length using star count analysis obtain scale length 
values between 3.5 and 4.5 kpc which is consistent with our 
analysis here, those works which have employed population 
se paration based on a kine matic or chemical sele ction, such 
as iBensbv et al.l (|201ll ) and ICheng et al.l (|2012d) , estimate 
a much lower scale length of around 2kpc. Understanding 
this difference will no doubt significantly improve definition 
of the astrophysical nature of the thick disk. 

While an increased hii reproduces the observed value 
of Ratio2, there is still a large disagreement in Ratiol (com- 
pare black line (data) and blue line (model) in Figure [22)l . 
Subtracting from the data the contribution of the high ve- 
locity substructure determined from Figure (red line), we 
see that the match improves for Ratiol too. 

Thus, from this simple toy-model analysis we deduce a 
value of ~ 4.0 kpc for the thick disk scale length. We find 
that a smooth double-exponential is, however, not a com- 
plete description of the data and a substantial overdensity, 
compared to this model, is apparent in the star counts, lo- 
calised near RA~313°. This overdensity matches the frac- 
tion of stars contributed by the high- velocity sub-structure 
identified in the previous section. 

In the next section, we analyse and characterise the 
stars in this sub-structure. 



4.5 Characterising the stars in the thick disk 
substructure 

In Figure[T9]we have already shown that the most significant 
contribution of stars which have high velocity compared to 
the model fit is restricted to just two fields in WGN (span- 
ning around 4° in RA) , beyond which the uniform thick disk 
model well reproduces the line of sight velocity distributions. 

Figure shows the colour-magnitude diagram of these 
two fields in the region R313 where the contribution of the 
high velocity stars is most significant. We see that both the 
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Figure 22. Comparison of SDSS photometric star count data in 
the WGN fields and the simple star count model described in the 
text. In the region RA^ 250° the star count data (black) and the 
simple model (blue) are in tolerable agreement. This indicates a 
photometric thick disk model with hfi,hz = 4.1,1.0 kpc is an 
adequate description of the Galaxy. For RA^ 300° the SDSS star 
count data (black) are significantly in excess of the model (blue), 
indicating an additional structure in the Galaxy. The red curve 
illustrates the SDSS count data after subtraction of a model for 
the excess counts described in the text and in Figure 1191 This 
"thick disk plus over-density" model is an adequate description 
of the underlying "normal" thick disk population. 
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Figure 23. Left : Colour-magnitude diagrams from WGN data 
for two fields, covering the region 312° < RA < 316° and 
-1° < Dec < 1°. Binsizcfor this CMD is Ag-r, Ag= 0.015,0.15. 
Right : CMD from SDSS DR,7 data for a larger region around the 
same fields covering 310° < RA < 318° and -3° < Dec < 3°. 
Binsize for this CMD is Ag - r,Ag= 0.01,0.1. The stars from 
the overdensity overlap with the thick disk turn off stars and are 
indistinguishable from thick disk. The colours denote low to high 
density from red to blue 



WGN and SDSS DR7 CMD look thick disk-like with the 
high-velocity stars overlapping with thick disk stars. No spe- 
cial feature in the CMD or extra turn-off is apparent. The 
substructure has a spread in apparent magnitudes as indi- 
cated by Figure 1241 but is preferentially seen at brighter 
magnitudes, being most prominent at magnitudes brighter 
than g ~ 17.3. The null hypothesis that the two curves 
(black and red) in Figure [24] are drawn from the same dis- 
tribution is rejected at the 5 percent level based on the Kol- 
mogorov Smirnov Test. The same comparison was done for 
RA~133 and the null hypothesis was accepted. 

Figure [24] also reveals that the substructure has metal- 
Ucity somewhat below the thick disk, with mean metallicity 
around -1.0 dex (black line). 

Is it reasonable for us to assume the excess stars are 



Figure 24. Unweighted line of sight velocity histograms are 
shown for all stars in three fields from WGN covering 312° < 
RA < 318°. The sample is divided into brighter and fainter 
than go = 17.3 The left panel shows the line of sight veloci- 
ties (histograms) with the blue dot-dash curve showing the ex- 
pected model distribution. The right panel shows the metallicity 
distributions. The excess stars are most apparent with velocities 
somewhat higher than that of the thick disk, with metallicity near 
[Fe/H]=-1, and at brighter magnitudes. 



in thick disk-like orbits? We looked at the Vr,V^,Vz dis- 
tribution of stars from this line of sight by cross matching 
the WGN data with proper motions from Stripe 82 based 
on position, being aware from our analysis above (Figure 
I16p that such comparisons are indicative rather than con- 
clusive. We used new high-prec ision values of SDSS star 
proper motions re-calibrated by IKoposov etahl l|2013h us- 
ing background galaxies. 613 matches were obtained. Stars 
with proper motion errors greater than 30mas/yr were re- 
moved. Figure [25] shows the velocity distributions for the 
region of interest (312° < RA < 317°) compared to distri- 
butions from nearby fields with no excess (RA > 317°). The 
and Vz distributions of the two sets of fields are indis- 
tinguishable. This is consistent with the excess stars being 
in thick disk-like orbits. Whereas, the mean Vr of the fields 
in R313 is offset towards positive values. 



4.6 Complementary and pre-discovery studies 

An independent source of kinematics in the region of in- 
terest comes from other SDSS studies, the SDSS spectra 
themselves. The SDSS DR7 spectra, selected applying the 
WGN selection function, hint at a cold high velocity struc- 
ture which is spatially and kinematically consistent with 
R313. Figure [26] shows line of sight velocities from SDSS 
DR7 for a region 312° < RA < 326° and -10° < Dec < 10°. 
Each panel is a 4° block in RA, the mid value of RA being 
shown at the top of each panel. The colours represent low 
to high density from red to blue while the dashed line shows 
the position of thick disk rotation (i.e V^=172 km/s) at a 
distance of 2 kpc for each direction and serves as a refer- 
ence. If velocities are consistent with being thick disk, we 
would expect the peak in density in each bin of declination 
to coincide with the dashed line. We can see that this is 
the case in most bins, apart from an excess of stars at rel- 
atively high Ygsr velocities (immediately to the right of 
the dashed line) in the first two panels at Dec~ 0°. This 
is consistent with the high velocity stars seen in the WGN 
data. For clarity, recall that these stars have observed radial 
velocity VheUo near zero (cf fig [3]). Any selection effects in 
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Figure 25. Galactocentric velocity distributions for the region with the star-count excess, R313 (black histogram), is compared to distri- 
butions from fields with RA > 317° which do not show a star count excess (red). The and Vz distributions are comparable, indicating 
a thick-disk-like velocity pattern, while the mean Vr of the fields in R313 is offset towards positive values. 
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Figure 26. Radial velocities from SDSS DR7 spectra for a region 312° < RA < 326° and -10° < Dec < 10°, applying the WGN selection 
function. Each panel is a 4° block in RA, the mid value of RA being shown at the top of each panel. The colours represent low to high 
density from red to blue. The dashed line shows the position of thick disk rotation (i.e. V0 = 172 km/s) at a distance of 2 kpc for each 
direction. We can see that there is an excess of stars at high Vcsij, that is here, close to zero Vhelio velocities (immediately to the right 
of the dashed line) in the first two panels, consistent with the excess stars seen in the WGN data. 



SDSS spectroscopy are unlikely to vary enough from field to 
field, within the tight cuts of the WGN selection function, 
to change this result. 

The overdensity observed by us is in a direction of the 
sky that is extremely busy, with previously identified sub- 
structures such as the Hercules Thick Disk Overdensity at 
1-2 kpc, Hercules Aquila Cloud at 10-20kpc and more recent 
sub-structure observed in RR Lyrae stars in Stripe 82 at 10- 
25kpc. 

The Hercules Thick Disk Overdensity, fir s t iden tified 
and studied in detail bv lLarsen fc Humphrevsl ('1996), and 
furthe r investigated in subs eq uent papers ([L arscn ct al.. 
(|2010[ ). iLarsen et all (l201ll ). lHumphrevs"et all (|201ll )) 
thereafter, is an excess of faint blue stars in the region 
I = 2Q° - 55° and 16| = 20° - 45°, thus overlapping with 
R313. These stars were estimated to have distances between 
1-2 kpc from the sun, consi stent with distanc e s esti mated 
by our method for our stars. [Humphreys et al] |20l3) addi- 
tionally also provide velocity measurements for their fields. 
They found that stars in Quadrant I (where the excess stars 
were found) had a Vlsr. of -28.6 ± 1.6 km/s above the plane 



and a Vlsr of -35.6 ± 2.0 km/s below the plane. They also 
state that going beyond 4 kpc along the line of sight, where 
the overdensity coincides with the density contours of the 
Galactic bar, the rotational rates become larger than stan- 
dard thick disk rotation. These kinematic parameters are 
in agreement with those we find here, though our sample 
he re is very much more sp atially concentrated than is that 
of lHumphrevs et al.1 (l201ll ). 

The Hercules Aquila Cloud jBelokurov et al.l [ioorl ) is 
thought to be a large structure further out in this same 
direction, centred on I ~ 40° and extending both above and 
below the plane up to 6 ~ 50°. The cloud was estimated 
to be more metal rich than -2.2 dex and have distances in 
the range of 10-20kpc. Their paper quotes a corresponding 
kinematic signature of Vgsr ~ 180km/s associated with 
this cloud which matches the kinematic signal we find in 
the WGN data. 

We als o looked at the density of RR Lyrae stars in 
Stripe 82 (jSesar et al.1 l2010l ). and towards the Hercules- 
Aquil a density enhancement. In Figure 11 of[s esar et al.l 
(|2010l ). the overdensity marked B is centred on similar direc- 
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tions as our overdensity, but extends between 5 and 25kpc. 
Cross-matching these stars with Sloan spectroscopy resulted 
in a handful of matches. While these stars showed a peak in 
[Fe/H] near -l.Sdex, we did not find a coherent velocity for 
stars in this structure. 

Given there are at least two known overdensities in this 
direction at distances beyond lOkpc, we investigated the 
possibility of the stars in our overdensity in fact being at 
distances beyond lOkpc, but being assigned grossly underes- 
timated distances by our by our analysis techniques. Look- 
ing at the predictions from the Besancjon model, we find 
that few sub-giants at large distances (further than lOkpc) 
pass through the WGN selection function and such stars 
fall in the very red, faint end of the WGN selection box. 
Our distance determination will indeed grossly underesti- 
mate the distances to these stars, putting them at 2-4 kpc. 
But this would require that all of the stars in the overden- 
sity be clearly concentrated in colour-magnitude space at 
the faint-red end. As seen from Figure [23] this is not the 
case. The overdensity stars instead look very thick-disk like. 



5 CONCLUSIONS 

We analyse the WGN medium resolution AAOmega spectro- 
scopic survey of some 10,000 stars, selected using SDSS pho- 
tometry to be main sequence turn-off stars in the thick disk- 
halo interface. The goal has been to investigate large-scale 
and small-scale kinematic, chemical and spatial properties 
of these two old stellar populations, and additionally also in- 
vestigate the region of transition/overlap between these two 
populations. 

Preparatory to the analysis, we developed a weighting 
method in colour-magnitude space for completeness correc- 
tion. Applying this correction, we obtained a completeness- 
corrected representation of the parent stellar population se- 
lected by the WGN selection function, thus preparing the 
data for analysis. Robust distance estimates were obtained 
by isochrone fitting, using either or both of SDSS gri and 
UKIDSS JHK photometry, for the stars, providing distances 
with a systematic uncertainty of less than 10%. 

We demonstrated how it is convenient to model the data 
in observable space directly. Propagated errors in derived pa- 
rameters, especially those in Vr,V^,Vz due to proper mo- 
tions, can be dominant, and hence deconvolution can de- 
grade the most accurate data, here radial velocities and 
[Fe/H] values, from any analysis. Hence, we analysed the 
kinematics of the thick disk using only observ able line of 
sight velocities, adopting the technique from ISirko et al.l 
(|2004l ) and limiting the use of derived parameters to that 
of distance alone. In this study, we obtain a global descrip- 
tion of the thick disk kinematics to demonstrate the scope 
of the data and the efficiency of the technique. The best 
fit thick disk kinematic parameters are determined to be 
V0 = 172 km/s and (j{(j),r,z) — (49,51,40) km/s which agrees 
reasonably well with values from previous studies. 

By comparing this global model description with the 
data, we identify a region near RA~313° which has higher 
line of sight Vgsr velocities than expected for the thick 
disk in this region. The same region is found also to have 
a local overdensity of stars in comparison to that expected 
for thick disk density models. The stars in the overdensity 



have a turn-off CMD similar to that of the thick disk and 
metallicities ~ — 1 dex, si milar to, but slight ly more metal 
poor than, the thick disk (|Carollo et al.|[2010l ). The stars in 
the excess are thus as old as are those of the thick disk. The 
overdensity decreases rapidly over a 4° increase in RA from 
312°. Its further extent is poorly determined as yet. 

This overdensity spatially may be associated in some 
way with the Hercules thic k disk overdensity. Line of sight 
velocity measurements from [Humphreys et al.l (|201ll ) match 
with the overdensity found in our WGN data, but they infer 
lower rotational velocities than thick disk for the stars in 
the overdensity while given our distance range we infer the 
opposite. We should have also expected the velocity signa- 
ture to continue till at least RA of 320° if associated with 
the Hercules thick disk overdensity which is a large struc- 
ture. This line of sight is towards the Hercules- Aquila halo 
cloud and a small spatial-scale overdensity in the RR Lyrae 
distribution. Obtaining reliable velocity estimates for these 
structures would be illuminating. 

We have shown that the method of modelling in observ- 
able space is very powerful not only to obtain a global de- 
scription of the Galactic stellar populations, but also to iden- 
tify and pick out substructure. We will extend this method 
in future work, to include scientifically important questions 
of spatial gradients in structure, kinematics and metallicity 
in the thick disk. 
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